{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The autoreload extension is already loaded. To reload it, use:\n",
      "  %reload_ext autoreload\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "'1.13.1'"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from __future__ import absolute_import, division, print_function\n",
    "%load_ext autoreload\n",
    "%autoreload 2\n",
    "import tensorflow as tf\n",
    "import numpy as np\n",
    "import experiments.MERA.modified_binary_mera as mod_bm\n",
    "import experiments.MERA.binary_mera as bm\n",
    "import experiments.MERA.binary_mera_lib as bml\n",
    "import tensornetwork as tn\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "#tf.enable_eager_execution()\n",
    "tf.enable_v2_behavior()\n",
    "tf.__version__"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# binary MERA\n",
    "A scale-invariant binary MERA is a tensor network of the following form:\n",
    "![title](MERA.png)\n",
    "\n",
    "The blue triangles are called isometries, and the green squares are called disentanglers.\n",
    "The microscopic degrees of freedom $i_n$ are at the bottom of the network.\n",
    "The following cell initializes a scale invariant binary MERA. \n",
    "Below `wC` is a list of `tf.Tensor` containing isometries in layers $\\tau$, \n",
    "and `uC` is a list of `tf.Tensor` with the disentanglers\n",
    "in layer $\\tau$. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The following cell performs a simple check: it calculates the expectatio value\n",
    "of the local Hamiltonian in each layer. The values should increase as $2^{\\tau}$.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "     Iteration: 1 of 10: E = -0.85743728, err = 0.4158022653439797 at D = 4 with 8 layers\n",
      "following numbers should all be 1/2\n",
      "[0.5 0.5 0.5 0.5 0.5 0.5 0.5 0.5]\n"
     ]
    }
   ],
   "source": [
    "def test_ascending_descending(chi=4,dtype = tf.float64):\n",
    "    import copy\n",
    "    wC, uC, rho_0 = bml.initialize_binary_MERA_identities(phys_dim=2, chi=4, dtype=dtype)\n",
    "    for n in range(5):\n",
    "        wC.append(copy.copy(wC[-1]))\n",
    "        uC.append(copy.copy(uC[-1]))\n",
    "    wC, uC,_,_ = bm.run_binary_mera_optimization_TFI(chis=[chi], niters=[10],opt_u_after=0 ,\n",
    "                                                       embeddings=['a'], dtype=dtype,\n",
    "                                                       wC=wC, uC=uC)\n",
    "    ham_0 = bml.initialize_TFI_hams(dtype)\n",
    "    rho = [0 for n in range(len(wC)+1)]\n",
    "    ham = [0 for n in range(len(wC)+1)]\n",
    "    rho[-1] = bml.steady_state_density_matrix(10, rho_0, wC[-1], uC[-1])\n",
    "    ham[0]= ham_0\n",
    "    print()\n",
    "    for p in range(len(rho)-2,-1,-1):\n",
    "        rho[p] = bml.descending_super_operator(rho[p+1], wC[p], uC[p])\n",
    "    for p in range(len(wC)):\n",
    "        ham[p+1] = bml.ascending_super_operator(ham[p], wC[p], uC[p])\n",
    "    energies = [tn.ncon([rho[p],ham[p]],[[1, 2, 3, 4, 5, 6],[1, 2, 3, 4, 5, 6]]) for p in range(len(rho))]\n",
    "    print('following numbers should all be 1/2')\n",
    "    print(np.array([energies[p]/energies[p+1] for p in range(len(energies) -1 )]))\n",
    "test_ascending_descending()    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## optimize \n",
    "The following cell runs an optimization to approximate the \n",
    "ground state of the transverse field Ising model with a MERA. \n",
    "The method returns list `uC` and `wC` of `tf.Tensor` objects, \n",
    "the disentanglers and the isometries, respectively."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "     Iteration: 391 of 400: E = -1.27323714, err = 0.0000024014674909 at D = 4 with 3 layerss"
     ]
    }
   ],
   "source": [
    "wC, uC,_,_ = bm.run_binary_mera_optimization_TFI(chis=[4,4], niters=[1000,400], \n",
    "                                                 embeddings=['a','a'], dtype=tf.float64, opt_u_after=40,\n",
    "                                                 opt_all_layers=[True, True], nsteps_steady_state=12)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "we can pass the tensors in `wC` and `uC` again to the routine to resume optimization"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "     Iteration: 291 of 300: E = -1.27323910, err = 0.0000004449449160 at D = 6 with 3 layers"
     ]
    }
   ],
   "source": [
    "wC, uC,_,_ = bm.run_binary_mera_optimization_TFI(chis=[6], niters=[300], \n",
    "                                                 embeddings=['p'], dtype=tf.float64,opt_u_after=0,\n",
    "                                                 opt_all_layers=[True], nsteps_steady_state=12,\n",
    "                                                 wC=wC, uC=uC)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Scaling dimensions\n",
    "From the optimized MERA, we can calculate scaling dimensions of the Ising CFT by diagonalizing\n",
    "the scale-invariant ascending super-operator:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEKCAYAAAD9xUlFAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAFsRJREFUeJzt3X+QXeV93/H3FyGBBHYUKophwVlS6HqonYwMg0mUetY4rQAT0NjgmPEvGDtKOijGqatUyrT2yJlMZCtDnNYUG9sYmTjYlKiyAoyVFHnjhAaKhKjFj26tgAAtGIwdISQW9MPf/nHPHi1itXsQOufs7n2/Zu7ce8499z7fZwT3s+ec5zwnMhNJkgCOarsASdLkYShIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpdHTbBbxWc+fOzTPOOKPtMhq1e/dujjvuuLbLaJR97g72uTmbNm16LjNPnGi7KRcKJ510Ehs3bmy7jEYNDAzQ39/fdhmNss/dwT43JyIer7Kdh48kSSVDQZJUMhQkSSVDQZJUMhQkSSVDQZJUMhQkSSVDQZJUiql2j+Z58+blkiVL2i6jUdu2baO3t7ftMhpln7uDfW7OihUrNmXmORNtN+VCoa+vLwcHB9suo1Fe9dkd7HN3aPGK5kqhMOWmuZCkqWjt5iFWrR9kaMcwPfdsYOnCPhbN72m7rFcxFCSpZms3D7F8zRaG9+4HYGjHMMvXbAGYdMHgiWZJqtmq9YNlIIwY3rufVesn36FwQ0GSavbUjuHXtL5NhoIk1eyUubNf0/o2GQqSusbazUMsWLmB05fdwYKVG1i7eaiRdpcu7GP2zBmvWDd75gyWLuxrpP3XorZQiIjTIuJ7EfFwRDwUEdeMsU1ExH+JiK0R8YOIeHtd9UjqbiMne4d2DJMcONnbRDAsmt/DH7/3bfQUewY9c2fzx+9926Q7yQz1jj7aB3wqM++PiDcAmyLibzLz4VHbXAicWTzeAVxfPEtdo62hiiPtPrVjmFPmzm50iGQbfR7vZG8T/V40v4dF83sm/bUZtYVCZj4NPF28fiEiHgF6gNGhcCnwjexcQXdPRMyNiJOLz0qNaesHsq2him0OkWyr7al0srdNjZxTiIheYD5w70Fv9QBPjlreXqyTGtPmYYW2hiq2OUSyrban0sneNtU+zUVEHA/8LfBHmbnmoPduB1Zm5t8Xy3cB/zEzNx603WJgMcCcOXPOvvzyy2utebLZs2cPs2bNaruMRjXZ53+YdxEvH33cq9Yfs283v/LcnbW2PXDSZRDx6jcy6X/mtmnXbpttP3PsaQy+8Rx+dtSBAyRH/WwffTs3ctJLT47zySOrrf+fV69e3f7cRxExE7gdWJ+Z147x/peBgcy8pVgeBPrHO3zk3Efdock+n77sDsb6vyCAx1a+p9a2F6zcwNAYhy965s7m7mXnT7t22267zfMoIyb73Ed1jj4K4GvAI2MFQmEd8JFiFNJ5wPOeT1DT2jys0NZQxTaHSLbZ9qL5Pdy97HweW/ke7l52/qQc/dO2OkcfLQA+DGyJiAeKdX8AvBkgM78E3AlcBGwFXgSuqrEeaUxLF/a94sQnNPsjBRwYidPQX6+j2236r+a2+qxq6hx99Pd09sDH2yaBq+uqQaqizR/IkfbbGKo40m4bpsrwzG7kLKkS7f5ASpOJ01xIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkq1hUJE3BgRz0bEg4d4vz8ino+IB4rHp+uqRZJUzdETbRARxwDvA3pHb5+Zn53gozcBXwS+Mc42f5eZF09YpSSpEROGAvAd4HlgE/By1S/OzO9HRO/hlSVJakNk5vgbRDyYmW89rC/vhMLtY30+IvqBvwS2A08B/yEzHzrE9ywGFgOceOKJZ996662HU86UtWvXLo4//vi2y2iUfe4O9rk573rXuzZl5jkTbVdlT+F/RcTbMnPLEahrtPuBX8jMXRFxEbAWOHOsDTPzBuAGgL6+vuzv7z/CpUxuAwMD2Ofpzz53h8ne5yonmn8N2BQRgxHxg4jYEhE/eL0NZ+bOzNxVvL4TmBkR817v90qSDl+VPYUL62g4It4EPJOZGRHn0gmon9TRliSpmglDITMfj4hfBv51servMvP/TPS5iLgF6AfmRcR24DPAzOI7vwRcBvy7iNgHDAMfyIlOcEiSalVlSOo1wG8Ba4pVfx4RN2Tmfx3vc5l5xQTvf5HOkFVJ0iRR5fDRx4B3ZOZugIj4HPAPwLihIEmaeqqcaA5g/6jl/cU6SdI0U2VP4evAvRHxP4rlRcDX6itJktSWKiear42IATpDUwGuyszNtVYlSWrFIUMhIt6YmTsj4gRgW/EYee+EzPxp/eVJkpo03p7CXwAX05nzaPRQ0SiWf7HGuiRJLThkKIzMXpqZpzdXjiSpTROOPoqIBRFxXPH6QxFxbUS8uf7SJElNqzIk9XrgxeKq5k8B/wjcXGtVkqRWVAmFfcX0E5cCX8zM64A31FuWJKkNVa5TeCEilgMfAt4ZEUdRzGEkSZpequwp/CadO659LDN/BJwKrKq1KklSK6pcvPYj4NpRy08w/n2XJUlTVJXRR++NiB9GxPMRsTMiXoiInU0UJ0lqVpVzCp8HfiMzH6m7GElSu6qcU3jGQJCk7lBlT2FjRHwbWEvnhDMAmbnm0B+RJE1FVULhjcCLwL8dtS45cCc2SdI0UWX00VVNFCJJal+V0Uf/MiLuiogHi+Vfioj/VH9pkqSmVTnR/BVgObAXIDN/AHygzqIkSe2oEgpzMvN/H7RuXx3FSJLaVSUUnouIf0Fxo52IuAx4utaqJEmtqDL66GrgBuAtETEEPEZncjxJ0jRTZfTRo8CvFzfaOSozX6i/LElSGyYMhYiYC3wE6AWOjggAMvMTtVYmSWpclcNHdwL3AFuAn9VbjiaLtZuHWLV+kKd2DHPK3NksXdjHovk9jbU7tGOYnns2NNaupI4qoXBsZv772ivRpLF28xDL12xheO9+AIZ2DLN8zRaAWn+g22pX0gFVRh/dHBG/FREnR8QJI4/aK1NrVq0fLH+YRwzv3c+q9YPTsl1JB0Tn9svjbBBxNfBHwA6KYalAZuYv1lzbmObNm5dLlixpo+nWbNu2jd7e3sbau+mlc4AY453kymM3Trt2J4um/50nA/vcnBUrVmzKzHMm2q5KKDwKnJuZzx2p4l6Pvr6+HBzsrr8cBwYG6O/vb6y9BSs3MLRj+FXre+bO5u5l50+7dieLpv+dJwP73JyIqBQKVQ4fbaUzS6q6xNKFfcyeOeMV62bPnMHShX3Tsl1JB1Q50bwbeCAivscr76fgkNRpauSkbtOjj0a3O7RjmJ4GRz1J6qgSCmuLh7rIovk9rfwYj7TbjYcVpMmgyhXNq5soRJLUvkOGQkTcmpnvj4gtHBh1VMrMX6q1MklS48bbU7imeL64iUIkSe07ZChk5tPF8+OH88URcSOdQHk2M986xvsB/BlwEZ3RTVdm5v2H05Yk6cg45JDUiHghInYe6lHhu28CLhjn/QuBM4vHYuD611J4N1i7eYgFKzdw5Xd3s2DlBtZuHmq7JEnT3Hh7Cm8AiIg/pHNTnZvpXG76QeDkib44M78fEb3jbHIp8I3sXD13T0TMjYiTR/ZQup3zAElqQ5WL1y7JzP+WmS9k5s7MvJ7OD/rr1QM8OWp5e7FOOA+QpHZUungtIj4IfIvOKKQr6FzQ1piIWEznEBNz5szhyiuvbLL5VgyddBnEq+cBGvqnF7ui/3v27OGmm25qu4xG2efuMNn7XGXuo146J4QX0AmFu4FPZua2Cb+889nbD3Gi+cvAQGbeUiwPAv0THT7qlrmPnAeo+y5es8/dYcrPfZSZ2zLz0sycl5knZuaiKoFQwTrgI9FxHvC85xMOcB4gSW2ocvjosETELUA/MC8itgOfAWYCZOaX6NzR7SIOTLh3VV21TEXOAySpDbWFQmZeMcH7CVxdV/vTgfMASWpaldFHkqQuMeGeQkSMdX/m54FNmfnAkS9JktSWKnsK5wC/Q+cagh7gt+lcqfyViPj9GmuTJDWsyjmFU4G3Z+YugIj4DHAH8E5gE/D5+sqTJDWpyp7CP2fUHdeAvcBJmTl80HpJ0hRXZU/hm8C9EfGdYvk3gL+IiOOAh2urTJLUuCp3XvvDiPgu8KvFqt/JzI3F6w/WVpkkqXFVr1O4Hxga2T4i3pyZT9RWlSSpFVWGpP4unauRnwH205k+OwFvxylJ00yVPYVrgL7M/EndxUiS2lVl9NGTdC5WkyRNc1X2FB4FBiLiDkYNQc3Ma2urSpLUiiqh8ETxmFU8JEnTVJUhqSuaKESS1L5DhkJEfCEzPxkRf0VntNErZOYltVYmSWrceHsKNxfPf9JEIZKk9h0yFDJzU/H8t82VI0lq03iHj7YwxmGjEZnpxWuSNM2Md/jo4saqkCRNCuMdPnq8yUIkSe2b8IrmiDgvIu6LiF0RsSci9kfEziaKkyQ1q8o0F18ErgB+CMwGPg5cV2dRkqR2VAkFMnMrMCMz92fm1+nco1mSNM1UmebixYiYBTwQEZ8HnqZimEiSppYqP+4fLrZbAuwGTgPeV2dRkqR2VNlTeA7Yk5kvASsiYgZwTL1lSZLaUGVP4S5gzqjl2cD/rKccSVKbqoTCsZm5a2SheD1nnO0lSVNUlVDYHRFvH1mIiLOB4fpKkiS1pco5hU8C/z0ingICeBPwm7VWJUlqRZWb7NwXEW8B+opVg5m5t96yJEltqDLNxeV0zis8CCwCvj36cJIkafqock7hP2fmCxHxa8C7ga8B19dbliSpDVVCYX/x/B7gK5l5BzCrvpIkSW2pEgpDEfFlOieX74yIYyp+TpI0xVT5cX8/sB5YmJk7gBOApbVWJUlqRZXRRy8Ca0YtP01nUjxJ0jTjYSBJUqnWUIiICyJiMCK2RsSyMd6/MiJ+HBEPFI+P11mPJGl8Va5oPizFbKrXAf8G2A7cFxHrMvPhgzb9dmYuqasOSVJ1de4pnAtszcxHM3MP8C3g0hrbkyS9TrXtKQA9wJOjlrcD7xhju/dFxDuB/wf8XmY+efAGEbEYWAxw4oknMjAwcOSrncR27dpln7uAfe4Ok73PdYZCFX8F3JKZL0fEbwOrgfMP3igzbwBuAOjr68v+/v5Gi2zbwMAA9nn6s8/dYbL3uc7DR0N0bt054tRiXSkzf5KZLxeLXwXOrrEeSdIE6gyF+4AzI+L0iJgFfABYN3qDiDh51OIlwCM11iNJmkBth48yc19ELKFzNfQM4MbMfCgiPgtszMx1wCci4hJgH/BT4Mq66pEkTazWcwqZeSdw50HrPj3q9XJgeZ01SJKq84pmSVLJUJAklQwFSVLJUJAklQwFSVLJUJAklQwFSVLJUJAklQwFSVLJUJAklQwFSVLJUJAklQwFSVLJUJAklQwFSVLJUJAklQwFSVLJUJAklQwFSVLJUJAklQwFSVLJUJAklQwFSVLJUJAklQwFSVLJUJAklY5uu4CpYO3mIVatH+SpHcOcMnc2Sxf2sWh+T9tlSdIRZyhMYO3mIZav2cLw3v0ADO0YZvmaLQAGg6Rpx8NHE1i1frAMhBHDe/ezav1gSxVJUn0iM9uu4TWZN29eLlmypLH2bnrpHCDGeCe58tiNjdSwbds2ent7G2lrsrDP3cE+N2fFihWbMvOcibabcqHQ19eXg4PN/ZW+YOUGhnYMv2p9z9zZ3L3s/EZqGBgYoL+/v5G2Jgv73B3sc3MiolIoePhoAksX9jF75oxXrJs9cwZLF/a1VJEk1ccTzRMYOZns6CNJ3cBQqGDR/B5DQFJX8PCRJKlkKEiSSoaCJKlUayhExAURMRgRWyNi2RjvHxMR3y7evzcieuusR5I0vtpCISJmANcBFwJnAVdExFkHbfYx4J8y8wzgT4HPTfS923b+jAUrN7B289CRLlmSul6dewrnAlsz89HM3AN8C7j0oG0uBVYXr28D3h0RY10+/Aoj8w8ZDJJ0ZNUZCj3Ak6OWtxfrxtwmM/cBzwP/rMqXO/+QJB15tU1zERGXARdk5seL5Q8D78jMJaO2ebDYZnux/I/FNs8d9F2LgcUAs950xtknf/QLnTcy6X/mtlrqn0z27NnDrFmz2i6jUfa5O9jn5qxevbrSNBd1Xrw2BJw2avnUYt1Y22yPiKOBnwN+cvAXZeYNwA0Ax5x8ZpliPT8/h5s+d9ORrXoScn6Y7mCfu0NbfV69evXEG1Hv4aP7gDMj4vSImAV8AFh30DbrgI8Wry8DNmTFXRfnH5KkI6+2PYXM3BcRS4D1wAzgxsx8KCI+C2zMzHXA14CbI2Ir8FM6wTGhHucfkqRa1Dr3UWbeCdx50LpPj3r9EnD5a/nO3jce1diU1ZLUbbyiWZJUMhQkSSVDQZJUMhQkSSVDQZJUqu2K5rpExAtAt81vMQ94bsKtphf73B3sc3N+ITNPnGijqXg7zsEql2pPJxGx0T5Pf/a5O0z2Pnv4SJJUMhQkSaWpGAo3tF1AC+xzd7DP3WFS93nKnWiWJNVnKu4pSJJqMqVCISIuiIjBiNgaEcvarqduEXFaRHwvIh6OiIci4pq2a2pCRMyIiM0RcXvbtTQhIuZGxG0R8X8j4pGI+JW2a6pbRPxe8d/0gxFxS0Qc23ZNR1pE3BgRzxY3ExtZd0JE/E1E/LB4/vk2axzLlAmFiJgBXAdcCJwFXBERZ7VbVe32AZ/KzLOA84Cru6DPANcAj7RdRIP+DPhuZr4F+GWmed8jogf4BHBOZr6VztT6labNn2JuAi44aN0y4K7MPBO4q1ieVKZMKADnAlsz89HM3AN8C7i05ZpqlZlPZ+b9xesX6PxYTOubSETEqcB7gK+2XUsTIuLngHfSubcImbknM3e0W1UjjgZmF3dcnAM81XI9R1xmfp/OfWJGuxQYuQXaamBRo0VVMJVCoQd4ctTydqb5D+RoEdELzAfubbeS2n0B+H3gZ20X0pDTgR8DXy8OmX01Io5ru6g6ZeYQ8CfAE8DTwPOZ+dftVtWYkzLz6eL1j4CT2ixmLFMpFLpWRBwP/CXwyczc2XY9dYmIi4FnM3NT27U06Gjg7cD1mTkf2M0kPKRwJBXH0S+lE4inAMdFxIfarap5xa2HJ93wz6kUCkPAaaOWTy3WTWsRMZNOIHwzM9e0XU/NFgCXRMQ2OocHz4+IP2+3pNptB7Zn5sge4G10QmI6+3Xgscz8cWbuBdYAv9pyTU15JiJOBiien225nleZSqFwH3BmRJweEbPonJha13JNtYqIoHOs+ZHMvLbteuqWmcsz89TM7KXz77shM6f1X5CZ+SPgyYjoK1a9G3i4xZKa8ARwXkTMKf4bfzfT/OT6KOuAjxavPwp8p8VaxjRlJsTLzH0RsQRYT2e0wo2Z+VDLZdVtAfBhYEtEPFCs+4Pi3teaPn4X+Gbxx86jwFUt11OrzLw3Im4D7qczwm4zk/wq38MREbcA/cC8iNgOfAZYCdwaER8DHgfe316FY/OKZklSaSodPpIk1cxQkCSVDAVJUslQkCSVDAVJUslQkCSVDAVJUslQkF6niOgt7oPwleIeAX8dEbPbrks6HIaCdGScCVyXmf8K2AG8r+V6pMNiKEhHxmOZOTIVySagt8VapMNmKEhHxsujXu9nCs0rJo1mKEiSSoaCJKnkLKmSpJJ7CpKkkqEgSSoZCpKkkqEgSSoZCpKkkqEgSSoZCpKkkqEgSSr9f4xDdtXNc6cfAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[-3.844111804577896e-15, 0.12507386444017216, 1.0002481919901671, 1.1242691320452483, 1.1249211992654748, 1.9863762926504265, 1.9975705618126158, 1.9991431763784262, 1.9994197140312209, 2.022228782482589, 2.0694707985718535]\n"
     ]
    }
   ],
   "source": [
    "k=11\n",
    "scdims = sorted(bml.get_scaling_dimensions(wC[-1], uC[-1], k=k))\n",
    "\n",
    "    \n",
    "\n",
    "%matplotlib inline\n",
    "plt.plot(range(len(scdims)),scdims,'o')\n",
    "plt.grid()\n",
    "plt.hlines(0.0,-2,12,linewidth=0.5)\n",
    "plt.hlines(0.125,-2,12,linewidth=0.5)\n",
    "plt.hlines(1.0,-2,12,linewidth=0.5)\n",
    "plt.hlines(1.125,-2,12,linewidth=0.5)\n",
    "plt.hlines(2,-2,12,linewidth=0.5)\n",
    "plt.hlines(2.125,-2,12,linewidth=0.5)\n",
    "plt.xlim([0,k])\n",
    "plt.xlabel('n')\n",
    "plt.ylabel('scaling dimension')\n",
    "plt.draw()\n",
    "plt.show()\n",
    "print(scdims)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# modified binary MERA\n",
    "the modified binary MERA is a simplified version of the binary MERA.  See webpage below:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 60,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<iframe src=https://www.tensors.net/p-mera width=950 height=1000></iframe>"
      ],
      "text/plain": [
       "<IPython.core.display.HTML object>"
      ]
     },
     "execution_count": 60,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import IPython\n",
    "url = 'https://www.tensors.net/p-mera'\n",
    "iframe = '<iframe src=' + url + ' width=950 height=1000></iframe>'\n",
    "IPython.display.HTML(iframe)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## optimize modified binary MERA"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "WARNING:tensorflow:From /Users/martin/PY3_TF_1.13/lib/python3.6/site-packages/tensorflow/python/ops/resource_variable_ops.py:642: colocate_with (from tensorflow.python.framework.ops) is deprecated and will be removed in a future version.\n",
      "Instructions for updating:\n",
      "Colocations handled automatically by placer.\n",
      "Iteration: 991 of 1000: E = -1.27322090, err = 0.0000186485689615 at D = 6 with 2 layers"
     ]
    }
   ],
   "source": [
    "#run a modified binary mera optimization for the TFI model\n",
    "energies, walltimes, wC, vV, uC = mod_bm.run_mod_binary_mera_optimization_TFI(chis=[4,6], niters=[1000, 1000], \n",
    "                                                                              embeddings=['p','a'],\n",
    "                                                                              dtype=tf.float64, verbose=1, \n",
    "                                                                              nsteps_steady_state=8,refsym=True)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY4AAAEOCAYAAACetPCkAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAH9tJREFUeJzt3Xl8HXW9//HXJ/seSpO0pVsKKUspQtu0BSqyiFr9UUBc2NcKVFn0/vTBBZer/n73Cup148cilQJlEQQRLAqi8hAqFWlToGwFaYOlBUr3Nkmz5/P7IycQc5tmTnLmzMk57+fjcR7JzJmZ84mM592Z7zLm7oiIiASVFXUBIiIyvCg4REQkLgoOERGJi4JDRETiouAQEZG4KDhERCQuCg4REYmLgkNEROKi4BARkbgoOEREJC45URcQhoqKCq+uro66DBGRYWXlypVb3L1yoO3SMjiqq6upq6uLugwRkWHFzNYF2U63qkREJC4KDhERiYuCQ0RE4qLgEBGRuCg4REQkLmkVHGY2z8wW7ty5M+pSRETSVloFh7s/4u6XlJeXR12KiEjaSqvgEBGR8Ck4REQkLgoOERGJi4JDRETiouAQEZG4KDhERCQuCg4REYmLgkNEROKi4BARkbikfHCY2f5mtsjMfh11LSIiEnJwmNltZrbJzF7us36umb1uZmvM7Oq9HcPd6919fph1iohIcGE/OvYO4Abgzp4VZpYN3Ah8DNgArDCzJUA2cG2f/S9y900h1ygiInEINTjcfamZVfdZPQtY4+71AGZ2H3CKu18LnBRmPSIiMnRRtHGMBdb3Wt4QW7dHZjbSzH4OTDOza/ay3SVmVmdmdZs3b05ctSIi8i/CvlU1ZO6+FVgQYLuFwEKA2tpaD7suEZFMFcUVx9vA+F7L42LrhkwPchIRCV8UwbECmGxmk8wsDzgDWJKIA+tBTiIi4Qu7O+69wDPAQWa2wczmu3sHcDnwOLAauN/dXwmzDhERSZywe1Wd2c/6R4FHE/15ZjYPmFdTU5PoQ4uISEzKjxyPh25ViYiEL62CQ0REwpdWwaFeVSIi4Uur4NCtKhGR8KVVcIiISPjSKjh0q0pEJHxpFRy6VSUiEr60Cg4REQmfgkNEROKS8rPjxqNn5Ph+1TU8/srGQR1jclUJ+1eWJLYwEZE0Yu7pNwN5/pjJPub8nw5u35wsFl80iyP3H5ngqkREUpuZrXT32gG3S8fgOPTwaX7/Y0/GvV97p/O1B1axcWcL911yJFPHqpFdRDJHRgdHbW2t19XVDWrfd3c289mbn6G5vZP7Lz2KmirdthKRzBA0ONQ43seY8kLu+cJsssw4d9GzvL2jOeqSRERSSloFR6IGAFZXFHPnRbNoau3g3FufZUtja4IqFBEZ/tIqOBI5AHDKfmXcfuFM3tnZzHmLlrOrpT0BFYqIDH9pFRyJNmPivtxybi1vbGpg/h0raG7rjLokEZHIKTgGcOyBlfz09GnUrdvOF+9ZSVtHV9QliYhESsERwP/60Bi+9+nDePL1zXz1gVV0dqVfTzQRkaDSauR4mM6cNYGdze1c99hrlBXk8J+nTsXMoi5LRCTp0io4eqYcqampCeX4C449gJ3N7dz85FrKC3O5au7BoXyOiEgqS6tbVcmYVv2qTxzEWbMncNOTa7nlqbWhfY6ISKpKqyuOZDAz/u8pU9nV3M61j71GeWEuZ8yaEHVZIiJJo+AYhOws48efP4LG1g6+/tBLlBXm8qnDxkRdlohIUqTVrapkysvJ4uazZzB9wgi+fN/zLP3H5qhLEhFJCgXHEBTmZbPogplMrirl0rtWsnLdtqhLEhEJnYJjiMoLc1l80SxGlxdw4e0rWP3urqhLEhEJlYIjASpL87lr/iyK83M4d9Fy/rmlKeqSRERCo+BIkHEjirhr/my63Dln0bNs3NkSdUkiIqFIq+BI1LTqg1VTVcLiC2exY3c75y56lu1NbZHUISISprQKjmQMABzIYePKufX8WtZt280Fty+nsbUjslpERMKQVsGRKo7cfyQ3nTWdl9/ZxcWL62hp13TsIpI+FBwhOXHKKP77cx/imfqtXHHv83R0ajp2EUkPCo4QfXraOL578qH86dX3uOrBF+nSdOwikgY05UjIzj+6mp3N7fz4T/+grCCXb8+bounYRWRYU3AkwRUn1LCzuZ1FT7/JPkW5fOXEA6MuSURk0BQcSWBmfONTh7CzuZ2f/vkNygtzuXDOpKjLEhEZFAVHkmRlGdeddhgNLe1895FXKSvI5TMzxkVdlohI3NQ4nkQ52Vn87IxpzKkZyVUPvsgfX9kYdUkiInFTcCRZQW42C8+t5bCx5Vx+7/P8be2WqEsSEYnLsAgOMzvVzH5hZr8ys49HXc9QFefncMeFM6keWcTFi+tYtX5H1CWJiAQWenCY2W1mtsnMXu6zfq6ZvW5ma8zs6r0dw90fdveLgQXA6WHWmyz7FOVx1/zZ7FuSx/m3L+eN9xqiLklEJJBkXHHcAcztvcLMsoEbgU8CU4AzzWyKmR1mZr/r86rqtes3Y/ulhVFlBdw9fza52Vmcs+hZ1m/bHXVJIiIDCj043H0p0PfReLOANe5e7+5twH3AKe7+kruf1Oe1ybp9H3jM3Z8Lu+ZkmjiymLvnz6alvYtzFz3LpgZNxy4iqW3A4DCzH5nZoQn+3LHA+l7LG2Lr+nMFcCLwWTNbsKcNzOwSM6szs7rNm4fX878PGl3K7RfO5L1drZy3aDk7d7dHXZKISL+CXHGsBhaa2bNmtsDMkj5nubtf7+4z3H2Bu/+8n20Wunutu9dWVlYmu8Qhmz5hBAvPm0H95iYuWryC3W2ajl1EUtOAweHut7r7HOA8oBp40cx+aWbHD+Fz3wbG91oeF1s3JFE/yGmojplcyfVnHsHzb21nwd3P0dahGXVFJPUEauOINWYfHHttAVYB/9vM7hvk564AJpvZJDPLA84AlgzyWO9LhQc5DdXcqWO47rQPsfQfm/m3X71Ap2bUFZEUM+CUI2b2E2Ae8ATwPXdfHnvr+2b2eoD97wWOAyrMbAPwbXdfZGaXA48D2cBt7v7KIP+GtPP5mePZ1dLOf/5+NaUFOVx72mGaUVdEUkaQuapeBL7p7k17eG/WQDu7+5n9rH8UeDTA5wdmZvOAeTU1NYk8bCS+cMz+7Njdzg1/WUN5US7XfPKQqEsSEQGCBccq4KA+/+LdCaxz95RqTHD3R4BHamtrL466lkT46scPZGdzO7c8VU95YS5fOm74B6KIDH9BguMmYDrdVx4GTAVeAcrN7Ivu/scQ64tLOl1xQPd07N89+VB2tbTzgz+8TnlhLmfPnhh1WSKS4YI0jr8DTIt1dZ0BTAPqgY8BPwizuHilQ+N4X1lZxn9/7nBOOLiKbz78Mo+seifqkkQkwwUJjgN7N1y7+6vAwe5eH15Z0ltudhY3nT2dmRP35asPrGLluu1RlyQiGSxIcLxqZjeb2bGx102xdfmAhjgnSUFuNrecO4Mx5QVcelcdb+9ojrokEclQQYLjfGAN8JXYqx64gO7QGMogwIQb7gMABzKiOI9F59fS2t7FFxbX0dSq0eUiknzm3v8As9jAvzvd/ezklTR0tbW1XldXF3UZoXnqH5u58PblfPSQUdxyzgyysjTGQ0SGzsxWunvtQNvt9YrD3TuBibHR3ZIijj2wkm+dNIU/vfoeP/zjgGMwRUQSKkh33HpgmZktAd4fBOjuPw6tKhnQBUdX88amRm5+ci2Tq0o4bfq4qEsSkQwRJDjWxl5ZQGm45QxNuo3j2JueMR7/3NLE1Q++xMSRRcyYuG/UZYlIBthrG8e/bGhW5O7D4hF16d7G0duO3W2ceuMyGlo6+O3lcxg3oijqkkRkmEpIG0fsQEeZ2avAa7Hlw2NdciUF7FOUx63nz6Sts7unVaN6WolIyIJ0x/0p8AlgK4C7rwI+EmZREp+aqhJuPGs6b2xq5Cv3aSp2EQlXoOdxuPv6Pqs6Q6hFhuAjB1byHydN4c+r3+MHj78WdTkiksaCNI6vN7OjATezXODLdD9ONuVkUuP4npx31ETe2NTALU/VM7mqlM/OUE8rEUm8IFccC4DLgLF0P971iNhyyknHSQ7jYWZ8e96hzKkZyTW/eZEV/9wWdUkikoaCPHN8i7uf7e6j3L3K3c9x963JKE7il5udxY1nTWfciCIuvWsl67cNi45wIjKMBOlVVWlmXzezhWZ2W88rGcXJ4HT3tKqlQz2tRCQEQW5V/RYoB/4M/L7XS1LYAZUl3HT2DNZsbuTL9z6vnlYikjBBgqPI3f/d3e939wd7XqFXJkP24ckVfGfeFJ54bRPf/4N6WolIYgQJjt+Z2adCryQB0n1a9cE496hqzj1yIguX1nN/Xd9e1SIi8QsSHF+mOzxazGyXmTWY2a6wCxuMTO9V1Z//mDeFD9dU8I2HXmL5m+ppJSJDE6RXVam7Z7l7gbuXxZbLklGcJEZPT6vxI4pYcLd6WonI0ATpVWVmdo6ZfSu2PN7MZoVfmiRSeVHu+z2t5i9eQUOLnvorIoMT5FbVTcBRwFmx5UbgxtAqktDsX1nCzefMYO3mJq5UTysRGaQgwTHb3S8DWgDcfTugJwIOU3NqKvjuyYfyl9c3c+2jKTlzjIikuCBzVbXHnj3u0D0gEOgKtSoJ1TlHTuSN9xq49ek3mTyqhNNnToi6JBEZRoJccVwPPARUmdl/AU8D3wu1Kgndt06awjGTK/jmwy/z93rNICMiwQXpVXUPcBVwLfAucKq7PxB2YYOhcRzB5WRnccNZ0xm/bxFfvHslb21VTysRCSbwo2OHk0x6dOxQvbmliVNvXEZVaT4Pfuloygpyoy5JRCKSsEfHSnqbVFHMzedM580tTXzp7ud4Z0dz1CWJSIpTcAhHH1DB9z59GH+v38qxP/wL1/zmJQ0SFJF+BelVJRng8zPHc3TNSH7+1FruX7GB++vW8+lpY7ns+BomVRRHXZ6IpJAB2zjM7DTg+0AVYLGXp/K0I2rjGJqNO1u4ZelafvnsW7R3dnHy4ftx+Qk11FSVRl2aiIQoaBtHkOBYA8xz92EzWkzBkRibGlq49a9vctcz62jp6ORTU8dw+Qk1HDImZf/NICJDkMjgWObucxJWWRIoOBJrW1Mbi56uZ/Hf1tHY2sHHp4ziyo9OZupYzUIskk4SGRw/A0YDDwOtPevd/TdDLTIsCo5w7Njdxu3L/snty95kV0sHxx9UyRUfncz0CSOiLk1EEiCRwXH7Hla7u1802OLCpuAI166Wdu56Zh23/rWe7bvbOWZyBVecMJlZk/aNurS96upy2jq76Ohy2ju6aO/qoqPT6ejsWd+93N5rm44up8sdB9wdd3Cn17rYevquj22Lx9b33q7nve59wtb9qUn4nCQNCYty5FllSR5zp46JsIJwJSw4hiMFR3I0tXZw99/X8Yu/1rOlsY3Zk/bl+IOrKMnPobSg55VLSX4OJfk5lBXkUpyfTU723nuBuzvN7Z00tXayu62D3W3dP7uXY7+3dbK7taOf5U6a2jrY3drJ7vbun01tHbS0a4o1GbrfXjaHw8fvE3UZoRhycJjZVe7+AzP7f+wh5N39yqGXOTAzO4TupxBWAE+4+80D7aPgSK7mtk7uXf4WC5fWs3FXy4DbF+Vld4dJQXegtHd6r2DooLm9M65/vRbmZlOcn01RXg5FedkU5WVTnN/9e3FeDoWx5cLcbPJyssjNNnKzs8jJziI3y7p/9qzL6nnvg2UzI8vAzDAgywyz7s/u+d0s9jvElnt+j+3LB9tZr+2yeg4UsuR8SvI+yJL3F72vua2Tj/3kKc6cNYHvnHxo0j8/GYIGx97GcfT0ohr0N7CZ3QacBGxy96m91s8FfgZkA7e6+3X9HSPWm2uBmWUBdwIDBockV2FeNhd9eBIXzqmmub2TxpYOdrV00NjaQUNLO40tHTS0dtDQ0tH9e0t793ut3ct5OVmxL/wcivOyKXr/Sz+2Lj+bwp733l+OhUJuNllZyf8Skcx00OhS3tLg2P6Dw90fif1cPITj3wHcQPcXPgCxKdpvBD4GbABWmNkSukPk2j77X+Tum8zsZOCLwF1DqEVCZmaxf/XnUKUeu5KGivNyaGzpiLqMyIU6ctzdl5pZdZ/Vs4A17l4PYGb3Aae4+7V0X53s6ThLgCVm9nvgl+FVLCLSv5KCHDZs13xuUUw5MhZY32t5AzC7v43N7DjgNCAfeHQv210CXAIwYYIeTCQiiVean0Nja3vUZURur8ERu610pbv/JEn1/A/u/iTwZIDtFgILobtxPNyqRCQTlRToVhUMMDuuu3cCZyb4M98GxvdaHhdbJyKS0orzc2hq7Yy6jMgFmVZ9mZndYGbHmNn0ntcQPnMFMNnMJplZHnAGsGQIx3ufngAoImEqyc+hrbOL1o7MDo8gbRxHxH7+n17rHDhhoB3N7F7gOKDCzDYA33b3RWZ2OfA43T2pbnP3V+Kquh+xnmCP1NbWXpyI44mI9FZa0P2V2djSQX5JdsTVRGfA4HD34wd7cHff420ud3+UvTR0D5aZzQPm1dTUJPrQIiKU5MeCo7WDkSX5EVcTnQFvVZlZuZn92MzqYq8fmVlKTovq7o+4+yXl5SlZnogMc8Wx4GjI8AbyIG0ctwENwOdjr13AniY+FBFJa6Wx4GhqzezgCNLGcYC7f6bX8nfN7IWwChoK3aoSkTCVFHxwqyqTBbniaDazD/csmNkcICWHTupWlYiEqXcbRyYLcsWxALizV7vGduD88EoSEUlNJWrjAAYeOZ4FHOTuh5tZGYC770pKZYOgW1UiEibdquo20MjxLuCq2O+7Ujk0QLeqRCRchbnZZJkax4O0cfzZzL5mZuPNbN+eV+iViYikGDOjJD9Ht6oCbHN67OdlvdY5sH/iyxERSW0l+TkZf6sqSBvHOe6+LEn1DInaOEQkbJohN1gbxw1JqmXI1MYhImHTFUewNo4nzOwzZqYHO4tIxispyFVwBNjmUuABoM3MdplZg5mldO8qEZGw7FOYy7amtqjLiNSAweHupe6e5e657l4WWy5LRnEiIqlmTHkBG3e14J65DxoNMjuumdk5Zvat2PJ4M5sVfmnx04OcRCRso8oKaOvoYvvuzH32eJBbVTcBRwFnxZYbgRtDq2gI1DguImEbU14AwLs7U3LKvqQIEhyz3f0yoAXA3bcDeaFWJSKSokbFguO9XS0RVxKdIMHRbmbZdA/6w8wqga5QqxIRSVEfXHEoOPbmeuAhoMrM/gt4GvheqFWJiKSoypJ8sgw2ZnBwBHnm+D1mthL4KGDAqe6+OvTKRERSUE52FpWl+QqOgbj7a8BrIdcyZJpyRESSYXR5IRvVxpEe1KtKRJJhdFlmX3GkVXCIiCTDmPJCNY6LiEhwo8sLaGztoKElMwcBKjhEROLU0yU3U29XKThEROI0prwQyNyxHAoOEZE4Zfq0IwoOEZE4jSorwExXHCIiElBeThYVJZnbJTetgkPTqotIsowpL+AdBcfwpwGAIpIso8sK2Kg2DhERCWq/fQp5d4euOEREJKDR5QU0tHbQ2NoRdSlJp+AQERmEqtJ8ADZl4GSHCg4RkUEYVdY9lmNTQ2vElSSfgkNEZBB6rjgy8RGyCg4RkUGoil1xbNYVh4iIBFFWkEN+TpauOEREJBgzY1RZgdo4REQkuKrSfF1xpCozKzazOjM7KepaRER66IojBGZ2m5ltMrOX+6yfa2avm9kaM7s6wKH+Hbg/nCpFRAansjSfTbsyLzhyQj7+HcANwJ09K8wsG7gR+BiwAVhhZkuAbODaPvtfBBwOvAoUhFyriEhcRpV1P0J2d1sHRXlhf52mjlD/UndfambVfVbPAta4ez2Amd0HnOLu1wL/41aUmR0HFANTgGYze9Tdu8KsW0QkiA9Gj7dSXaHgCNNYYH2v5Q3A7P42dvdvAJjZBcCW/kLDzC4BLgGYMGFComoVEelXz+jx93a1UF1RHHE1yTMsGscB3P0Od//dXt5f6O617l5bWVmZzNJEJENVlcWuODKsgTyK4HgbGN9reVxs3ZDpQU4ikkyjSj+44sgkUQTHCmCymU0yszzgDGBJIg6sBzmJSDKVFeaQl5OVcdOOhN0d917gGeAgM9tgZvPdvQO4HHgcWA3c7+6vhFmHiEgYukePZ94gwLB7VZ3Zz/pHgUcT/XlmNg+YV1NTk+hDi4jsUVVp5g0CHDaN40HoVpWIJFsmXnGkVXCIiCSbrjiGOfWqEpFkqyrLp6Glg+a2zqhLSZq0Cg7dqhKRZOvpkrsxg25XpVVwiIgk2/h9iwB4a9vuiCtJnrQKDt2qEpFkqx7ZHRzrtjZFXEnypFVw6FaViCRbZWk+RXnZvLlFwSEiIgGYGdUji1m7WcEhIiIBHTKmjNXv7oq6jKRJq+BQG4eIROGQMaVsbmjNmDmr0io41MYhIlE4bGz3d84L63dEXElypFVwiIhE4fDx+5CXk8XyN7dGXUpSKDhERIaoIDeb2okjeGL1Jtw96nJCp+AQEUmAeYfvR/2WJl55J/0bydMqONQ4LiJR+eTU0eRmGw89n5AHmqa0tAoONY6LSFT2KcrjE4eO5lcr1rOtqS3qckKVVsEhIhKlr5w4md1tHVz/xBtRlxIqBYeISILUVJVy1uwJLH7mn/y9Pn17WCk4REQS6JpPHsLEfYu4/JfP8dbW9JwxV8EhIpJAxfk5LLpgJh1dznm3Pcu7O5ujLinh0io41KtKRFLBAZUlLDp/Jlsb2zjtpr/x2sb06qKbVsGhXlUikipmTBzBry49is4u55QblnHv8rfSZnBgWgWHiEgqmbJfGb+/8hhmVu/LNb95ifNuW0795saoyxoyBYeISIgqS/O586JZfGfeFF54awef+OlSvrPkFd7ZMXzbPixdLp16q62t9bq6uqjLEBH5F5sbWvnRH1/n1ys3YAanTRvHuUdNZOrY1Li9bmYr3b12wO0UHCIiybVh+25ueaqeB1aup6W9i8PGlvO52nHMPXQ0VWUFkdWl4FBwiEiK29nczm9feJtfPvsWr21swAxqJ45g7tQxHHtgJQdUFmNmSatHwaHgEJFhwt15Y1Mjj720kcdefpfXNjYAMKosnzkHVHDUASOZNmEf9q8oISsrvCBRcCg4RGSYemvrbpat3cKyNVv429qt70+aWJKfw6H7lXH4+H2YXFXCAVUlHFBRQnlRbkI+NyODw8zmAfNqamoufuON9J5kTEQyQ1eXs2ZzI6vW7+Clt3eyasNOVr+zi7bOrve3GVmcx6SKYsaNKORLx9dw4KjSQX1WRgZHD11xiEg6a+/sYv223dRvbqJ+SyP1m5t4c0sTb+9o5ufnzBh0L62gwZEzqKOLiEhkcrOz2L+yhP0rS4BRSf98DQAUEZG4KDhERCQuCg4REYmLgkNEROKi4BARkbgoOEREJC4KDhERiYuCQ0RE4pKWI8fNbDOwDigH+nsAeX/vVQBbQiotEfb2N6XCsQdzjKD7BNluoG36e39v+6XyORHm+ZCI44d5PgTZdijvZ+J3xER3rxxwK3dP2xewMN73gLqo6x7s35QKxx7MMYLuE2S7gbbZy3/3vZ0rKXtOhHk+JOL4YZ4PQ/nvHfC/u74j+nml+62qRwb5XioLs+5EHHswxwi6T5DtBtqmv/d1PoRz/DDPhyDbDuV9nRP9SMtbVUNhZnUeYJIvyRw6J6Q3nQ9qHN+ThVEXIClH54T0lvHng644REQkLrriEBGRuCg4REQkLgoOERGJi4JjAGZWbGaLzewXZnZ21PVItMxsfzNbZGa/jroWSQ1mdmrs++FXZvbxqOtJhowMDjO7zcw2mdnLfdbPNbPXzWyNmV0dW30a8Gt3vxg4OenFSujiOR/cvd7d50dTqSRLnOfEw7HvhwXA6VHUm2wZGRzAHcDc3ivMLBu4EfgkMAU408ymAOOA9bHNOpNYoyTPHQQ/HyQz3EH858Q3Y++nvYwMDndfCmzrs3oWsCb2L8o24D7gFGAD3eEBGfq/V7qL83yQDBDPOWHdvg885u7PJbvWKOiL8ANj+eDKAroDYyzwG+AzZnYzw3cKAonfHs8HMxtpZj8HppnZNdGUJhHp7zviCuBE4LNmtiCKwpItJ+oCUp27NwEXRl2HpAZ330r3vWwRANz9euD6qOtIJl1xfOBtYHyv5XGxdZKZdD5IXzonYhQcH1gBTDazSWaWB5wBLIm4JomOzgfpS+dETEYGh5ndCzwDHGRmG8xsvrt3AJcDjwOrgfvd/ZUo65Tk0Pkgfemc2DtNcigiInHJyCsOEREZPAWHiIjERcEhIiJxUXCIiEhcFBwiIhIXBYeIiMRFwSGyF2b2t9jPajM7K8HH/vqePksk1Wkch0gAZnYc8DV3PymOfXJig8b6e7/R3UsSUZ9IMumKQ2QvzKwx9ut1wDFm9oKZ/ZuZZZvZD81shZm9aGaXxrY/zsz+amZLgFdj6x42s5Vm9oqZXRJbdx1QGDvePb0/KzZN9w/N7GUze8nMTu917CfN7Ndm9pqZ3WNmltz/RUQ0O65IUFfT64ojFgA73X2mmeUDy8zsj7FtpwNT3f3N2PJF7r7NzAqBFWb2oLtfbWaXu/sRe/is04AjgMOBitg+S2PvTQMOBd4BlgFzgKcT/+eK9E9XHCKD83HgPDN7AXgWGAlMjr23vFdoAFxpZquAv9M9u+pk9u7DwL3u3unu7wFPATN7HXuDu3cBLwDVCflrROKgKw6RwTHgCnd//F9WdreFNPVZPhE4yt13m9mTQMEQPre11++d6P/DEgFdcYgE0wCU9lp+HPiimeUCmNmBZla8h/3Kge2x0DgYOLLXe+09+/fxV+D0WDtKJfARYHlC/gqRBNC/VkSCeRHojN1yugP4Gd23iZ6LNVBvBk7dw35/ABaY2WrgdbpvV/VYCLxoZs+5+9m91j8EHAWsAhy4yt03xoJHJHLqjisiInHRrSoREYmLgkNEROKi4BARkbgoOEREJC4KDhERiYuCQ0RE4qLgEBGRuCg4REQkLv8fsM7vBJZFA5UAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.loglog(np.array(energies)+4/np.pi)\n",
    "plt.xlabel('iteration')\n",
    "plt.ylabel('error in energy')\n",
    "plt.draw()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
